Fourier methods for the perturbed harmonic oscillator in linear and nonlinear Schr\"odinger equations

We consider the numerical integration of the Gross-Pitaevskii equation with a potential trap given by a time-dependent harmonic potential or a small perturbation thereof. Splitting methods are frequently used with Fourier techniques since the system can be split into the kinetic and remaining part, and each part can be solved efficiently using Fast Fourier Transforms. To split the system into the quantum harmonic oscillator problem and the remaining part allows to get higher accuracies in many cases, but it requires to change between Hermite basis functions and the coordinate space, and this is not efficient for time-dependent frequencies or strong nonlinearities. We show how to build new methods which combine the advantages of using Fourier methods while solving the timedependent harmonic oscillator exactly (or with a high accuracy by using a Magnus integrator and an appropriate decomposition).Comment: 12 pages of RevTex4-1, 8 figures; substantially revised and extended versio

High order structure preserving explicit methods for solving linear-quadratic optimal control problems

[EN] We consider the numerical integration of linear-quadratic optimal control problems. This problem requires the solution of a boundary value problem: a non-autonomous matrix Riccati differential equation (RDE) with final conditions coupled with the state vector equation with initial conditions. The RDE has positive definite matrix solution and to numerically preserve this qualitative property we propose first to integrate this equation backward in time with a sufficiently accurate scheme. Then, this problem turns into an initial value problem, and we analyse splitting and Magnus integrators for the forward time integration which preserve the positive definite matrix solutions for the RDE. Duplicating the time as two new coordinates and using appropriate splitting methods, high order methods preserving the desired property can be obtained. The schemes make sequential computations and do not require the storrage of intermediate results, so the storage requirements are minimal. The proposed methods are also adapted for solving linear-quadratic N-player differential games. The performance of the splitting methods can be considerably improved if the system is a perturbation of an exactly solvable problem and the system is properly split. Some numerical examples illustrate the performance of the proposed methods.The author wishes to thank the University of California San Diego for its hospitality where part of this work was done. He also acknowledges the support of the Ministerio de Ciencia e Innovacion (Spain) under the coordinated project MTM2010-18246-C03. The author also acknowledges the suggestions by the referees to improve the presentation of this work.Blanes Zamora, S. (2015). High order structure preserving explicit methods for solving linear-quadratic optimal control problems. Numerical Algorithms. 69:271-290. https://doi.org/10.1007/s11075-014-9894-0S27129069Abou-Kandil, H., Freiling, G., Ionescy, V., Jank, G.: Matrix Riccati equations in control and systems theory. Basel, Burkhäuser Verlag (2003)Al-Mohy, A.H., Higham, N.J.: Computing the Action of the Matrix Exponential, with an Application to Exponential Integrators. SIAM. J. Sci. Comp. 33, 488–511 (2011)Anderson, B.D.O., Moore, J.B.: Optimal control: linear quadratic methods. Dover, New York (1990)Ascher, U.M., Mattheij, R.M., Russell, R.D.: Numerical solutions of boundary value problems for ordinary differential equations. Prentice-Hall, Englewood Cliffs (1988)Bader, P., Blanes, S., Ponsoda, E.: Structure preserving integrators for solving linear quadratic optimal control problems with applications to describe the flight of a quadrotor. J. Comput. Appl. Math. 262, 223–233 (2014)Basar, T., Olsder, G.J.: Dynamic non cooperative game theory, 2nd Ed, SIAM, Philadelphhia (1999)Blanes, S., Casas, F.: On the necessity of negative coefficients for operator splitting schemes of order higher than two. Appl. Num. Math. 54, 23–37 (2005)Blanes, S., Casas, F., Farrés, A., Laskar, J., Makazaga, J., Murua, A.: New families of symplectic splitting methods for numerical integration in dynamical astronomy. Appl. Numer. Math. 68, 58–72 (2013)Blanes, S., Casas, F., Oteo, J.A., Ros, J.: The Magnus expansion and some of its applications. Phys. Rep. 470, 151–238 (2009)Blanes, S., Casas, F., Ros, J.: High order optimized geometric integrators for linear differential equations. BIT 42, 262–284 (2002)Blanes, S., Diele, F., Marangi, C., Ragni, S.: Splitting and composition methods for explicit time dependence in separable dynamical systems. J. Comput. Appl. Math. 235, 646–659 (2010)Blanes, S., Moan, P.C.: Practical symplectic partitioned Runge-Kutta and Runge-Kutta-Nystrm methods. J. Comput. Appl. Math. 142, 313–330 (2002)Blanes, S., Ponsoda, E.: Magnus integrators for solving linear-quadratic differential games. J. Comput. Appl. Math. 236, 3394–3408 (2012)Brif, C., Chakrabarti, R., Rabitz, H.: Control of quantum phenomena: past, present and future. New J. Phys. 12, 075008(68pp) (2010)Cruz, J.B., Chen, C.I.: Series Nash solution of two person non zero sum linear quadratic games. J. Optim. Theory Appl. 7, 240–257 (1971)Dieci, L., Eirola, T.: Positive definitness in the numerical solution of Riccati differential quations. Numer. Math. 67, 303–313 (1994)Engwerda, J.: LQ dynamic optimization and differential games. Wiley (2005)Hairer, E., Lubich, C., Wanner, G.: Geometric Numerical Integration. Structure-Preserving Algorithms for Ordinary Differential Equations (2nd edition). Springer Series in Computational Mathematics, 31. Springer-Verlag (2006)Hochbruck, M., Ostermann, A.: Exponential integrators. Acta Numerica 19, 209–286 (2010)Horn, R.A., Johnson, C.R.: Matrix Analysis. Cambridge University Press, New York (1985)Iserles, A., Munthe-Kaas, H.Z., Nørsett, S.P., Zanna, A.: Lie group methods. Acta Numerica 9, 215–365 (2000)Iserles, A., Nørsett, S.P.: On the solution of linear differential equations in Lie groups. Phil. Trans. R. Soc. Lond. A 357, 983–1019 (1999)Jódar, L., Ponsoda, E.: Non-autonomous Riccati-type matrix differential equations: existence interval, construction of continuous numerical solutions and error bounds. IMA. J. Num. Anal. 15, 61–74 (1995)Jódar, L., Ponsoda, E., Company, R.: Solutions of coupled Riccati equations arising in differential games. Control. Cybern. 24, 117–128 (1995)Kaitala, V, Pohjola, M. In: Carraro, Filar (eds.) : Sustainable international agreement on greenhouse warming. A game theory study. Control and Game Theoretic Models of the Environment, pp 67–87. Birkhauser, Boston (1995)Keller, H.B.: Numerical solution of two point boundary value problems. In: CBMS-NSF Regional Conference Series in Applied Mathematics, Vol. 24. SIAM, Philadelphia (1976)McLachlan, R.I.: Composition methods in the presence of small parameters. BIT 35, 258–268 (1995)McLachlan, R.I., Quispel, R.: Splitting Methods. Acta Numer. 11, 341–434 (2002)Moler, C.B., Van Loan, C.F.: Nineteen Dubious Ways to Compute the Exponential of a Matrix, twenty-five years later. SIAM Rev. 45, 3–49 (2003)Na, T.Y.: Computational methods in engineering boundary value problems. In: Mathematics in Science and Engineering, Vol. 145. Accademic Press, New York (1979)Palao, J.P., Kosloff, R.: Quantum computing by an optimal control algorithm for unitry transformations. Phys. Rev. Lett. 28 (2002)Peirce, A.P., Dahleh, M.A., Rabitz, H.: Optimal control of quantum-mechanical systems: existence, numerical approximation, and applications. Phys. Rev. A 37, 4950–4967 (1988)Reid, W.T.: Riccati Differential Equations. Academic, New York (1972)Sanz-Serna, J.M., Calvo, M.P.: Numerical Hamiltonian Problems. Chapman & Hall, London (1994)Sidje, R.B.: Expokit: a software package for computing matrix exponentials. ACM Trans. Math. Software 24, 130–156 (1998)Speyer, J.L., Jacobson, D.H.: Primer on optimal control theory. SIAM, Philadelphia (2010)Starr, A.W., Ho, Y.C.: Non-zero sum differential games. J. Optim. Theory and Appl 3, 179–197 (1969)Zhu, W., Rabitz, H.: A rapid monotonically convergent iteration algorithm for quantum optimal control ever the expectation value of a positive definite operator. J. Chem. Phys. 109, 385–391 (1998

Carbon Dioxide Production in Animal Houses: A literature review

This article deals with carbon dioxide production from farm animals; more specifically, it addresses the possibilities of using the measured carbon dioxide concentration in animal houses as basis for estimation of ventilation flow (as the ventilation flow is a key parameter of aerial emissions from animal houses). The investigations include measurements in respiration chambers and in animal houses, mainly for growing pigs and broilers. Over the last decade a fixed carbon dioxide production of 185 litres per hour per heat production unit, hpu (i.e. 1000 W of the total animal heat production at 20 oC) has often been used. The article shows that the carbon dioxide production per hpu increases with increasing respiration quotient. As the respiration quotient increases with body mass for growing animals, the carbon dioxide production per heat production unit also increases with increased body mass. The carbon dioxide production is e.g. less than 185 litres per hour per hpu for weaners and broilers and higher for growing finishing pigs and cows. The analyses show that the measured carbon dioxide production is higher in full scale animal houses than measured in respiration chambers, due to differences in manure handling. In respiration chambers there is none or very limited carbon dioxide contribution from manure; unlike in animal houses, where a certain carbon dioxide contribution from manure handling may be foreseen. Therefore, it is necessary to make a correction of data from respiration chambers, when used in full scale animal buildings as basis for estimation of ventilation flow. Based on the data reviewed in this study, we recommend adding 10% carbon dioxide production to the laboratory based carbon dioxide production for animal houses with slatted or solid floors, provided that indoor manure cellars are emptied regularly in a four weeks interval. Due to a high and variable carbon dioxide production in deep straw litter houses and houses with indoor storage of manure longer than four weeks, we do not recommend to calculate the ventilation flow based on the carbon dioxide concentration for these houses

Simulation of bimolecular reactions: Numerical challenges with the graph Laplacian

[EN] An important framework for modelling and simulation of chemical reactions is a Markov process sometimes known as a master equation. Explicit solutions of master equations are rare; in general the explicit solution of the governing master equation for a bimolecular reaction remains an open question. We show that a solution is possible in special cases. One method of solution is diagonalization. The crucial class of matrices that describe this family of models are non-symmetric graph Laplacians. We illustrate how standard numerical algorithms for finding eigenvalues fail for the non-symmetric graph Laplacians that arise in master equations for models of chemical kinetics. We propose a novel way to explore the pseudospectra of the non-symmetric graph Laplacians that arise in this class of applications, and illustrate our proposal by Monte Carlo. Finally, we apply the Magnus expansion, which provides a method of simulation when rates change in time. Again the graph Laplacian structure presents some unique issues: standard numerical methods of more than second-order fail to preserve positivity. We therefore propose a method that achieves fourth-order accuracy, and maintain positivity.We thank the organisers and delegates of the Canberra 2019 EMAC conference for helpful discussions about graph Laplacians. SM thanks the Australian Research Council Centre of Excellence fot Mathematical ans Statistical Frontiers (ACEMS). The work of SB was funded by Ministerio de Economía, Industria y Competitividad (Spain) through project MTM2016-77660-P (AEI/FEDER, UE).Macnamara, S.; Blanes Zamora, S.; Iserles, A. (2020). Simulation of bimolecular reactions: Numerical challenges with the graph Laplacian. The ANZIAM Journal. 61:1-16. https://doi.org/10.21914/anziamj.v61i0.15169S1166

High-order splitting methods for separable non-autonomous parabolic equations

We consider the numerical integration of non-autonomous separable parabolic equations using high order splitting methods with complex coefficients (methods with real coeffi- cients of order greater than two necessarily have negative coefficients). We propose to consider a class of methods that allows us to evaluate all time-dependent operators at real values of the time, leading to schemes which are stable and simple to implement. If the system can be considered as the perturbation of an exactly solvable problem and the flow of the dominant part is advanced using real coefficients, it is possible to build highly efficient methods for these problems. We show the performance of this class of methods on several numerical examples and present some new improved schemesThe authors thank the referees for their suggestions to improve the presentation of this work. The work of Sergio Blanes has been supported by Ministerio de Ciencia e Innovacion (Spain) under project MTM2010-18246-C03 and the Ministerio de Educacion, Cultura y Deporte, under Programa Nacional de Movilidad de Recursos Humanos del Plan Nacional de I-D+i 2008-2011 (PRX12/00547). The work of Muaz Seydaoglu has been supported by the Turkish Council of High Education through a grant for visiting the Instituto de Matematica Multidisciplinar at the Polytechnic University of Valencia where this work was carried out.Seydaoglu, M.; Blanes Zamora, S. (2014). High-order splitting methods for separable non-autonomous parabolic equations. Applied Numerical Mathematics. 84:22-32. https://doi.org/10.1016/j.apnum.2014.05.004S22328

Solving the Schrödinger eigenvalue problem by the imaginary time propagation technique using splitting methods with complex coefficients

The Schrodinger eigenvalue problem is solved with the imaginary time propagation technique. The separability of the Hamiltonian makes the problem suitable for the application of splitting methods. High order fractional time steps of order greater than two necessarily have negative steps and cannot be used for this class of diffusive problems. However, there exist methods which use fractional complex time steps with positive real parts which can be used with only a moderate increase in the computational cost. We analyze the performance of this class of schemes and propose new methods which outperform the existing ones in most cases. On the other hand, if the gradient of the potential is available, methods up to fourth order with real and positive coefficients exist. We also explore this case and propose new methods as well as sixth-order methods with complex coefficients. In particular, highly optimized sixth-order schemes for near integrable systems using positive real part complex coefficients with and without modified potentials are presented. A time-stepping variable order algorithm is proposed and numerical results show the enhanced efficiency of the new methods.We wish to acknowledge Ander Murua and Joseba Makazaga for providing the methods T869 and V869. This work has been partially supported by Ministerio de Ciencia e Innovacion (Spain) under Project MTM2010-18246-C03 and by a grant from the Qatar National Research Fund (NPRP) #NPRP 5-674-1-114. P.B. also acknowledges the support through the FPU fellowship AP2009-1892.Bader, P.; Blanes Zamora, S.; Casas, F. (2013). Solving the Schrödinger eigenvalue problem by the imaginary time propagation technique using splitting methods with complex coefficients. Journal of Chemical Physics. 139(12):124117-124117. https://doi.org/10.1063/1.4821126S12411712411713912Chin, S. A., Janecek, S., & Krotscheck, E. (2009). Any order imaginary time propagation method for solving the Schrödinger equation. Chemical Physics Letters, 470(4-6), 342-346. doi:10.1016/j.cplett.2009.01.068Janecek, S., & Krotscheck, E. (2008). A fast and simple program for solving local Schrödinger equations in two and three dimensions. Computer Physics Communications, 178(11), 835-842. doi:10.1016/j.cpc.2008.01.035Roy, A. K., Gupta, N., & Deb, B. M. (2001). Time-dependent quantum-mechanical calculation of ground and excited states of anharmonic and double-well oscillators. Physical Review A, 65(1). doi:10.1103/physreva.65.012109Auer, J., Krotscheck, E., & Chin, S. A. (2001). A fourth-order real-space algorithm for solving local Schrödinger equations. The Journal of Chemical Physics, 115(15), 6841-6846. doi:10.1063/1.1404142Lehtovaara, L., Toivanen, J., & Eloranta, J. (2007). Solution of time-independent Schrödinger equation by the imaginary time propagation method. Journal of Computational Physics, 221(1), 148-157. doi:10.1016/j.jcp.2006.06.006Trefethen, L. N., & Bau, D. (1997). Numerical Linear Algebra. doi:10.1137/1.9780898719574Aichinger, M., & Krotscheck, E. (2005). A fast configuration space method for solving local Kohn–Sham equations. Computational Materials Science, 34(2), 188-212. doi:10.1016/j.commatsci.2004.11.002SHENG, Q. (1989). Solving Linear Partial Differential Equations by Exponential Splitting. IMA Journal of Numerical Analysis, 9(2), 199-212. doi:10.1093/imanum/9.2.199Suzuki, M. (1991). General theory of fractal path integrals with applications to many‐body theories and statistical physics. Journal of Mathematical Physics, 32(2), 400-407. doi:10.1063/1.529425Chin, S. A. (1997). Symplectic integrators from composite operator factorizations. Physics Letters A, 226(6), 344-348. doi:10.1016/s0375-9601(97)00003-0Omelyan, I. P., Mryglod, I. M., & Folk, R. (2002). Construction of high-order force-gradient algorithms for integration of motion in classical and quantum systems. Physical Review E, 66(2). doi:10.1103/physreve.66.026701Chin, S. A., & Chen, C. R. (2002). Gradient symplectic algorithms for solving the Schrödinger equation with time-dependent potentials. The Journal of Chemical Physics, 117(4), 1409-1415. doi:10.1063/1.1485725Goldman, D., & Kaper, T. J. (1996). Nth-Order Operator Splitting Schemes and Nonreversible Systems. SIAM Journal on Numerical Analysis, 33(1), 349-367. doi:10.1137/0733018Blanes, S., & Casas, F. (2005). On the necessity of negative coefficients for operator splitting schemes of order higher than two. Applied Numerical Mathematics, 54(1), 23-37. doi:10.1016/j.apnum.2004.10.005E. Hairer, C. Lubich, and G. Wanner, Geometric Numerical Integration. Structure-Preserving Algorithms for Ordinary Differential Equations, 2nd ed. (Springer, Berlin, 2006), p. 81.P.V. Koseleff, Ph.D. thesis, École Polytechnique, 1993.Chin, S. A. (2005). Structure of positive decompositions of exponential operators. Physical Review E, 71(1). doi:10.1103/physreve.71.016703Chin, S. A., & Krotscheck, E. (2005). Fourth-order algorithms for solving the imaginary-time Gross-Pitaevskii equation in a rotating anisotropic trap. Physical Review E, 72(3). doi:10.1103/physreve.72.036705Bader, P., & Blanes, S. (2011). Fourier methods for the perturbed harmonic oscillator in linear and nonlinear Schrödinger equations. Physical Review E, 83(4). doi:10.1103/physreve.83.046711Chambers, J. E. (2003). Symplectic Integrators with Complex Time Steps. The Astronomical Journal, 126(2), 1119-1126. doi:10.1086/376844Bandrauk, A. D., Dehghanian, E., & Lu, H. (2006). Complex integration steps in decomposition of quantum exponential evolution operators. Chemical Physics Letters, 419(4-6), 346-350. doi:10.1016/j.cplett.2005.12.006BANDRAUK, A. D., & LU, H. (2013). EXPONENTIAL PROPAGATORS (INTEGRATORS) FOR THE TIME-DEPENDENT SCHRÖDINGER EQUATION. Journal of Theoretical and Computational Chemistry, 12(06), 1340001. doi:10.1142/s0219633613400014Castella, F., Chartier, P., Descombes, S., & Vilmart, G. (2009). Splitting methods with complex times for parabolic equations. BIT Numerical Mathematics, 49(3), 487-508. doi:10.1007/s10543-009-0235-yHansen, E., & Ostermann, A. (2009). High order splitting methods for analytic semigroups exist. BIT Numerical Mathematics, 49(3), 527-542. doi:10.1007/s10543-009-0236-xBlanes, S., Casas, F., Chartier, P., & Murua, A. (2012). Optimized high-order splitting methods for some classes of parabolic equations. Mathematics of Computation, 82(283), 1559-1576. doi:10.1090/s0025-5718-2012-02657-3Blanes, S., Casas, F., & Ros, J. (2000). Celestial Mechanics and Dynamical Astronomy, 77(1), 17-36. doi:10.1023/a:1008311025472McLachlan, R. I. (1995). Composition methods in the presence of small parameters. BIT Numerical Mathematics, 35(2), 258-268. doi:10.1007/bf01737165Laskar, J., & Robutel, P. (2001). Celestial Mechanics and Dynamical Astronomy, 80(1), 39-62. doi:10.1023/a:1012098603882Sakkos, K., Casulleras, J., & Boronat, J. (2009). High order Chin actions in path integral Monte Carlo. The Journal of Chemical Physics, 130(20), 204109. doi:10.1063/1.3143522Blanes, S., Casas, F., & Ros, J. (1999). Celestial Mechanics and Dynamical Astronomy, 75(2), 149-161. doi:10.1023/a:1008364504014Blanes, S., Casas, F., Farrés, A., Laskar, J., Makazaga, J., & Murua, A. (2013). New families of symplectic splitting methods for numerical integration in dynamical astronomy. Applied Numerical Mathematics, 68, 58-72. doi:10.1016/j.apnum.2013.01.003Blanes, S., Casas, F., & Ros, J. (2001). High-order Runge–Kutta–Nyström geometric methods with processing. Applied Numerical Mathematics, 39(3-4), 245-259. doi:10.1016/s0168-9274(00)00035-

The Magnus expansion and some of its applications

Approximate resolution of linear systems of differential equations with varying coefficients is a recurrent problem shared by a number of scientific and engineering areas, ranging from Quantum Mechanics to Control Theory. When formulated in operator or matrix form, the Magnus expansion furnishes an elegant setting to built up approximate exponential representations of the solution of the system. It provides a power series expansion for the corresponding exponent and is sometimes referred to as Time-Dependent Exponential Perturbation Theory. Every Magnus approximant corresponds in Perturbation Theory to a partial re-summation of infinite terms with the important additional property of preserving at any order certain symmetries of the exact solution. The goal of this review is threefold. First, to collect a number of developments scattered through half a century of scientific literature on Magnus expansion. They concern the methods for the generation of terms in the expansion, estimates of the radius of convergence of the series, generalizations and related non-perturbative expansions. Second, to provide a bridge with its implementation as generator of especial purpose numerical integration methods, a field of intense activity during the last decade. Third, to illustrate with examples the kind of results one can expect from Magnus expansion in comparison with those from both perturbative schemes and standard numerical integrators. We buttress this issue with a revision of the wide range of physical applications found by Magnus expansion in the literature.Comment: Report on the Magnus expansion for differential equations and its applications to several physical problem

New families of symplectic splitting methods for numerical integration in dynamical astronomy

We present new splitting methods designed for the numerical integration of near-integrable Hamiltonian systems, and in particular for planetary N-body problems, when one is interested in very accurate results over a large time span. We derive in a systematic way an independent set of necessary and sufficient conditions to be satisfied by the coefficients of splitting methods to achieve a prescribed order of accuracy. Splitting methods satisfying such (generalized) order conditions are appropriate in particular for the numerical simulation of the Solar System described in Jacobi coordinates. We show that, when using Poincar\'e Heliocentric coordinates, the same order of accuracy may be obtained by imposing an additional polynomial equation on the coefficients of the splitting method. We construct several splitting methods appropriate for each of the two sets of coordinates by solving the corresponding systems of polynomial equations and finding the optimal solutions. The experiments reported here indicate that the efficiency of our new schemes is clearly superior to previous integrators when high accuracy is required.Comment: 24 pages, 2 figures. Revised version, accepted for publication in Applied Numerical Mathematic

New efficient numerical methods to describe the heat transfer in a solid medium

The analysis of heat conduction through a solid with heat generation leads to a linear matrix differential equation with separated boundary conditions. We present a symmetric second order exponential integrator for the numerical integration of this problem using the imbedding formulation. An algorithm to implement this explicit method in an efficient way with respect to the computational cost of the scheme is presented. This method can also be used for nonlinear boundary value problems if the quasilinearization technique is considered. Some numerical examples illustrate the performance of this method. © 2010 Elsevier Ltd.The authors acknowledge the support of the Generalitat Valenciana through the project GV/2009/032 and the Ministerio de Ciencia e Innovacion (Spain) under projects MTM2007-61572 and MTM2009-08587 (co-financed by the ERDF of the European Union).Ponsoda Miralles, E.; Blanes Zamora, S.; Bader, P. (2011). New efficient numerical methods to describe the heat transfer in a solid medium. Mathematical and Computer Modelling. 54(7-8):1858-1862. doi:10.1016/j.mcm.2010.11.067S18581862547-

Improved high order integrators based on the Magnus expansion

Abstract. We build high order efficient numerical integration methods for solving the linear differential equationẊ = A(t)X based on the Magnus expansion. These methods preserve qualitative geometric properties of the exact solution and involve the use of single integrals and fewer commutators than previously published schemes. Sixth-and eighth-order numerical algorithms with automatic step size control are constructed explicitly. The analysis is carried out by using the theory of free Lie algebras